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Spin dynamics in a doped-Mott-insulator superconductor 

W.Q. Chen and Z.Y. Weng 
Center for Advanced Study, Tsinghua University, Beijing 100084, China 

We present a systematic study of spin dynamics in a superconducting ground state, which itself 
is a doped-Mott-insulator and can correctly reduce to an antiferromagnetic (AF) state at half- 
filling with an AF long-range order (AFLRO). Such a doped Mott insulator is described by a 
mean-field theory based on the phase string formulation of the t — J model. We show that the 
well-known spin wave excitation in the AFLRO state at half-filling evolves into a resonancelike 
„Cj ■ peak at a finite energy in the superconducting state, which is located around the AF wave vectors. 

D ' The width of such a resonancelike peak in momentum space decides a spin correlation length scale 

which is inversely proportional to the square root of doping concentration, while the energy of the 
fvj ■ resonancelike peak scales linearly with the doping concentration at low doping. These properties 

r*0 ' are consistent with experimental observations in the high-T c cuprates. An important prediction of 

the theory is that, while the total spin sum rule is satisfied at different doping concentrations, the 
weight of the resonancelike peak does not vanish, but is continuously saturated to the weight of 
the AFLRO at zero-doping limit. Besides the low-energy resonancelike peak, we also show that the 
high-energy excitations still track the spin wave dispersion in momentum space, contributing to a 
significant portion of the total spin sum rule. The fluctuational effect beyond the mean-field theory 
is also examined, which is related to the broadening of the resonancelike peak in energy space. In 
Mh' particular, we discuss the incommensurability of the spin dynamics by pointing out that its visibility 

is strongly tied to the low-energy fluctuations below the resonancelike peak. We finally investigate 
the interlayer coupling effect on the spin dynamics as a function of doping, by considering a bilayer 
system. 
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I. INTRODUCTION 



The measurement of spin dynamics in the cuprate superconductors is uniquely important. This is because the spin 

degrees of freedom constitute the predominant part of the low- lying electronic degrees of freedom, i.e., 1 — S per 

site, as compared to the charge degrees of freedom at small hole concentration, 5 per site. Such a large imbalance 

£T J i between the spin and charge numbers are usually regarded as a key indicator that the underlying system is a doped 

>y-s Mott insulator 1 . On general grounds, the corresponding spin dynamics is expected to be distinctly different from 

fj . a conventional BCS superconductor. The latter is based on the Fermi-liquid description in which the elementary 

excitations are quasiparticles that carry both charge and spin. An extreme case is at half-filling, where the whole 

charge degrees of freedom get frozen at low energy and only the spin degrees of freedom remain intact in the cuprates, 

^^ ■ whose dynamics is well characterized by the Heisenberg model 2 . 

Experimentally, anomalous properties of spin dynamics have been observed throughout the cuprate family. The 

parent compound at half-filling is a Mott insulator in which spins form AFLRO below a Neel temperature TV- The 

elementary excitation is a gapless bosonic Goldstone mode, i.e., the spin wave in the ordered phase. AFLRO and the 

spin-wave excitation disappear beyond some critical concentration of holes introduced into the system. Except for 

some residual signature of spin waves at high energies, the low-lying spin-wave-type excitation is completely absent 

^j , once the system becomes a superconductor. It is replaced by a resonancelike peak at a doping-dependent energy 

around the AF wave vector Qaf = (7r,7r), as observed first in the optimally doped YBCO compound 3,4 , where 

. i-H ', the dynamic spin susceptibility function measured by inelastic neutron scattering shows a sharp peak at aj ros = 41 

j^ ■ meV, whose width is comparable to the resolution limit of the instruments. Similar resonancelike peak has also been 

observed in the underdoped YBCO compounds (where the resonancelike peak persists into the pseudogap phase 

above the superconducting transition), Tl-based 6 and Bi-based 7 compounds. In the LSCO compound, although no 

such a sharp peak has been found, the low-lying spin excitation is nonetheless non spin- wave-like, which may be still 

regarded as a very broad peak in energy space 8 . With much sharper linewidth in momentum space, doping-dependent 

incommensurate splittings around Qaf have been clearly identified in LSCO 9,10 . Similar incommensurability, even 

though not as prominent as in LSCO, has been also established in underdoped YBCO recently 11-13 . 

Theoretically, a great challenge is how to naturally connect the spin dynamic at half-filling with that in the super- 
conducting phase in which the doping concentration can be as low as a few percent per Cu site. That is, although the 
low-energy, long-wavelength behavior may change qualitatively in the superconducting phase, the number of spins in 
the background is still quite close to half-filling, which far exceeds the number of doped holes. Physically it is very 



hard to imagine that the short-range, high-energy spin correlations would be changed completely by a few percent 
to ten percent doping. However, in a BCS superconductor, the upper spin energy scale is usually set by the Fermi 
energy e/ 14 , such that in the local spin susceptibility one has to integrate over the frequency up to ef in order to 
recover the correct sum rule of 1 — S spin per site. Normally e/ is much larger than J. Thus, why there should be a 
gigantic increase in the upper spin energy in the doped case, compared to the half-filling, poses a serious challenge to 
any approach based on the d-wave BCS-typc theory in which the spin dynamic is solely contributed by quasiparticlc 
excitations. Experimentally the upper energy scale exhibited in the dynamic spin susceptibility is set by ~ 2J (J is 
the superexchangc coupling) at half-filling, in consistency with the prediction by the Hcisenbcrg model, and is slightly 
reduced in the optimal-doped superconducting phase 8 . No trace of any other new high-energy scale has been ever 
reported in the doped regime in spin channels. 

As for the low-energy feature, like the resonancelike peak structure observed in the experiments, theoretical proposals 
are ranged from the RPA fluctuations in the particle-hole channel within the framework of BCS 14 or generalized BCS 
theories 15,16 to some novel mechanism of the so-called it mode in the particle-particle channel in the SO(5) theory 17 , 
which is coupled to the particle-hole channel in the superconducting phase. An important question, not being properly 
addressed yet, is what is the connection, if any, of such a resonancclike spin mode with the spin wave in the zero- 
doping limit. Namely, how a few percent of doped holes can continuously reshape a spin- wave excitation into a 
non-propagating local mode, with an AFLRO turning into short-range spin correlations. This question and the 
previous high-energy one constitute two of most fundamental issues in an approach based on doped Mott insulators. 

In this paper, we put forward a systematic description of the evolution of spin dynamics as a function of doping 
in a doped-Mott-insulator superconductor. It is described by a bosonic resonating-valence-bond (RVB) mean-field 
theory 18 based on the phase-string formulation 19 of the t — J model. At half-filling, the mean-field theory reduces 
to the Schwinger-boson mean- field state 20 , which well characterizes AFLRO and spin- wave excitations in the ground 
state. At finite doping, the mean-field theory depicts how the spin dynamics is influenced by the doping effect in 
going into the superconducting state. In particular, we show how a resonancelike peak centered around Qaf emerges 
out of spin waves from the AFLRO phase. A unique prediction for experiment is that the weight of the resonancclike 
peak continuously evolves into that of the AFLRO in the zero-doping limit. On the other hand, the total weight of 
the dynamic susceptibility function, which extends slightly over ~ 2 J in energy, still satisfies the sum rule that the 
total spin number is 1 — 5 per site. 

In this unified mean-field description, doping-dependent rcsonancelike energy and spin correlation length are quan- 
titatively determined. Besides the low-energy resonancelike peak structure near Qaf> there still exists a high-energy 
spectrum whose envelope roughly tracks the spin wave dispersion as a residual effect in the superconducting phase. 
We also consider some leading fluctuational effect beyond the mean-field theory on the lineshape of the spectral 
function, and discuss the incommensurability and its visibility in this framework. We finally introduce the interlayer 
superexchangc coupling and investigate how the spin dynamics changes in the even and odd channels for a double-layer 
system. Comparisons with the experimental measurements, mostly by inelastic neutron scattering, are made. 

The remainder of the paper is organized as follows. In Sec. II, a systematic study of spin dynamics in the bosonic 
RVB mean-field state for the single-layer system is presented. In Sec. Ill, fluctuational effects beyond the mean-field 
theory, due to the charge density fluctuations, are discussed. In Sec. IV, the interlayer coupling for a bilayer system 
is considered. Finally, a summary is given in Sec. V. 

II. SPIN DYNAMICS IN MEAN FIELD DESCRIPTION 
A. Bosonic RVB state at half-filling 

Spin dynamics of the cuprates at half-filling is well described by the two-dimensional (2D) AF Hcisenberg model. 
Although a conventional spin-wave theory is quite successful in understanding the low-lying excitation spectrum of 
the Hcisenberg Hamiltonian, to make the theory applicable or modifiable to the cases without AFLRO, like at finite 
temperatures or in doped regimes, we shall use the Schwinger-boson formulation as our starting point at half-filling. 

The mean-field theory 20 based on the Schwinger-boson formulation can characterize the AFLRO and spin-wave 
excitation fairly well in the ground state. Its mean-field wavefunction under the Gutzwiller projection will have the 
same form 21 as the variational bosonic RVB wavefunctions proposed by Liang, Doucot, and Anderson 22 . The latter 
can produce very accurate variational energies as well as the AF magnetization for the Hcisenberg model, indicating 
that the state correctly captures both long-range and short-range spin correlations. Such an approach is thus called 
bosonic RVB description, which is to be generalized to finite doping in the next subsection. In the following, we 
briefly review some basic equations of the bosonic RVB mean-field theory at half-filling. 

In the Schwinger-boson formulation, the spin operators can be expressed by the Schwinger-boson operator bi„ as 
follows 



S+ = (-l)%b ii , (1) 



(note that a staggered sign factor (—1)* is explicitly introduced here in contrast to the original definition 20 ), and 
S~ = (S^Y, while Sf = ^ CT crb icr bi a . The Schwinger bosons satisfy the constraint ^ CT b icT bi a = 1. The mean-field 
state is characterized by the bosonic RVB order parameter 

Ag = £>„&,-_,), (2) 

IT 

which leads to the following effective Hamiltonian, obtained from the half-filling t — J (Hciscnbcrg) model: 

T A s / \ 

h s = —^ Y. b l b U + H - c - + const - + A E b l b * - N > ( 3 ) 

(ij)a \ ia / 

where the last term involves a Lagrangian multiplier A to enforce the global constraint of total spinon number, 
Ei„blb ia = N. 

The mean-field Heisenberg Hamiltonian (3) can be straightforwardly diagonalizcd by the Bogoliubov transformation 

bio = E w ko-(*)( u k7k<r - Wk7k-cr)' ( 4 ) 



as 



h s = J2 e ^L^- ( 5 ) 



ka 



Here, Wk CT (i) = ~~7f! e%a , and the coherent factors, u^ and i>k, are given by 



1 / A sgn(Ck) /A 



where £k = — JAq(cos k x a + cos fcj,a) and 



E* = \V-& (7) 



Finally, in a self-consistent manner, the RVB order parameter Ag and the Lagrangian multiplier A are determined by 
the following self-consistent equations 

|Ag| 2 = ^y^-coth^!i, (8) 

1 ol 2N ^ JE k 2 ' V ; 

k 

2 = IyA co th^ + n ' ( 9 ) 

N ^ E u 2 BC ' y ' 

k/0 K 

in which n b BC is the contribution from the Bose condensation of the Schwinger bosons, leading to an AFLRO, which 
happens if E^ becomes gapless. Note that (3 = 1/T and the AFLRO disappears (n BC = 0) at a finite temperature T. 

B. Bosonic RVB description at finite doping 

Although AF correlations at half-filling are well captured by the mean- field Hamiltonian H s in (3), the doping effect 
on the spin background is a highly nontrivial issue. 

Based on the phase-string formulation 19 , which is an exact reformulation by sorting out the most singular doping 
effect, i.e., the phase string effect in the t — J model, a generalized mean- field Hamiltonian describing the spin degrees 
of freedom can be obtained 18 as follows 



H * = ~^y~ E b \A-° eiaA% + H - c - + const - + A (E 6 U* -{i-S)n\. (io) 

(ij)cr \ icr / 

Compared to the half-filling case, H s in (10) differs from (3) mainly by the emergence of a gauge field A\- defined 
on the nearest-neighboring (NN) link (ij), satisfying the following constraint 

E 4 = * E »?> (") 

where c is a, say, counter-clockwise-oriented close loop and f2 c is the area enclosed by c. On the right-hand-side (rhs) 
of (11), nj 1 denotes the number operator of doped holes at site I. Therefore, the doping effect explicitly enters in (10) 
through the gauge field A\a as if each hole carries a fictitious 7r fluxoid as seen by spinons in H s . In (10), the bosonic 
RVB order parameter is given by 



A . = £/ e -^ ff6 \ (12) 



JVN 



for NN sites i and j. At half filling, because there is no hole, it is obvious that A\, = 0, and A s reduces back to Aq 
defined in (2). 

Note that the doping concentration 5 also enters H s through the Lagrangian multiplier A which implements the 
global condition J^. b\ a bi a = (1 — 5)N. But at low doping, the effect of missing spins represented by such a term will 
be far less dramatic than the topological gauge field A\- . The latter reflects the singular phase string effect 19 induced 
by the hopping of doped holes on the AF spin background. 

Corresponding to H s in (10), the spin operators in the phase string formulation 19 read 

St = (-\) i J*h\,b ii , (13) 

S^~ = (S i ; )t, and Sf = ^ CT crb ilJ bi ( y 1 respectively. Compared to the Schwinger-boson formulation in (1), an extra phase 
<J>f appears in (13), which satisfies <j>^ — ^ = 2A^ (ij £ spin sites) and ensures the spin rotational symmetry of the 
effective Hamiltonian (10). 

Equation (10) is by nature a gauge model. But in the superconducting ground state, due to the Bose condensation 
of bosonic holons in the bosonic RVB theory 18 , the spin Hamiltonian H s will become quite simplified as A\a can be 
approximately treated as describing a uniform flux with a strength 

E D 4 = ^- ( 14 ) 

Then we can introduce the following Bogoliubov transformation to diagonalize (10), just like (4) in diagonalizing 

(3), 

bio = ^W m «j(i)( u m7mu - Vm'fL-o-)- ( 15 ) 

m 

With a standard procedure, we obtain 

•ffs = E E ^Lalma, (16) 

ma 

with the spinon spectrum 

E m = v/A^e- (17) 

In this scheme, £ m and Lo ma (i) — uj* n _ <J (i) are eigenvalues and eigenfunctions of the following equation 

JA. S h 

£, m 0Jma{i) = ^~ JZ ^""WO')- ( 18 ) 

j=NN(i) 



and the coherent factors, u m and v m , are given by 



Finally, A and A s can be determined by the self-consistent equations 

1 ' 2NJ ^ E m 2 ' 

m m 



2 -^^E^-th: 



''BC- 



(19) 

(20) 
(21) 



Here n b BC is the contribution of the condensation of spinons, if an AFLRO exists like in the half-filling case. 

The above mean-field formulation is essentially the same as the one outlined in Ref. 18 . For simplicity and clarity, 
here we have not explicitly included an approximate doping-correction factor in (20) («1- 25) as we will be mainly 
concerned with the evolution of spin dynamics at low doping. Such additional corrections from doped holes can be 
always incorporated by simply replacing the supcrcxchange coupling J with a doping-dcpcndent J e ff which is quickly 
reduced at higher doping concentrations. A spin feedback effect from the hopping term is not included cither, which 
results in a shift of A to A m in E m ls , without qualitatively changing the physical consequences. 



C. Spin dynamics in superconducting ground state 
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FIG. 1: The density of states (DOS) of the mean-field spinon spectrum E m at doping 8 — 0.125. Inset: the DOS in the AF 
state at half filling. 



1. Excitation spectrum E„ 



According to the mean-field scheme outlined above, we can numerically determine the mean-field 'spinon' spectrum 
E m defined in (17). 



As an example, we solve the eigenequation (18) and self-consistent equations (20) and (21) at doping concentration 
S = 0.125. The chemical potential A is found to be 1.819J while the RVB order parameter A s is 0.993. In contrast, 
at half filling, the results are A = 2.316J and Aq = 1.158. 

In Fig. 1, the density of states (DOS) of the spectrum E m for 6 = 0.125 is shown in the main panel, while the 
half-filling case is plotted in the inset for comparison. The figure shows that the two spectra are qualitatively very 
different. At half-filling, the spectrum is continuous and gapless, with a large density of states at the maximal energy 
which is slightly above 2 J. In the superconducting state, the spectrum becomes discretized levels. This discrete levels 
are due to the fact that the spectrum £ m as the solution of (18) has a Hofstadter spectrum as the result of the vector 
potential Afj given in (14). Note that the distribution of the Landau-level-like structure in Fig. 1 remains uneven, 
which reflects the fact that the average density of states increases with energy, as seen at half-filling. The maximal 
energy is slightly less than 2 J at 5 = 0.125. 

It is important to note that there is a gap between the lowest discrete level and zero energy, which is ~ 0.265J for 
S = 0.125. There will no more spinon Bosc condensation n h BC ^ such that the AFLRO no longer exists. 

2. Dynamic spin susceptibility 



After diagonalizing the effective Hamiltonian H s , the spin susceptibility can be obtained straightforwardly. Due to 
the spin rotational invariance 18 , one may only consider the z-component susceptibility, which can be derived based 
on the Matsubara Green's function — (T t Sj(t)S£(0)). With the standard procedure outlined in Ref. 18 , the imaginary 
part of the dynamic spin susceptibility at zero temperature is given by 

X"(Q,«) - ^ E Crnm'iQ) ( - ~ t^ - - l) sgn(w)<J(M - E m - E m ,), (22) 

mm 

where 



Cw(Q) = ^ E e iQ,(xrIj) ^(i)4(j>^(i)^(j'). (23) 

The discrete energy levels of E m illustrated in Fig.l will show up in x"(Q- w). We plot the positions of these peaks 
in x"(Q, lu) in energy and momentum space, as well as the FWHM (full width of half maximum) in momentum space, 
in Fig. 2. The momentum scan in Fig. 2(a) is along (it, q) direction and is along the diagonal (q, q) in Fig. 2(b). One 
sees that each discrete energy corresponds to a finite width in momentum as depicted by a finite bar. 

For comparison, the spin-wave peak positions at half-filling are shown as dotted curves in Fig. 2. At S = 0.125, 
although the spin excitations are no longer propagating modes, as evidenced by the flat (dispcrsionlcss) small bars at 
discrete energies, the envelope of the overall spectrum at high energies still approximately track the dispersion of the 
spin wave at half-filling, with a slightly softened spin-wave velocity. Note that there are actually some more peaks at 
even higher energies than in Fig. 2, but their weight is much reduced due to the coherent factors in %"(Q,w) (see the 
local spin susceptibility below). 

Fig. 2 clearly depicts how the spin excitations in the superconducting state continuously evolves from the spin- 
wave picture at half-filling. The remnant high-energy spin wave signature at finite doping is a very unique feature 
in this approach. Recently, such a high-energy spin wave feature has been reported 2 in underdoped YBa2Cu30g.5 
compound. 

In the following, we turn our attention to the lowest peak in Fig. 2, which has the largest weight as marked by the 
darkest FWHM bar. 

3. Resonancelike peak around AF wave-vector Qaf 



Let us consider two special momenta, Qo = (0,0) and Qaf = (it, it). For Q = Qo, with the relation 
I] i w mC r(«) w m'o'(*) = fimm', wc have C TOTO '(Q ) = jjSmm' such that 

x"(Qo, w) = -^ J2 ( A2 ~f 2 ~^ ) s § n (^(M - 2E ™) = °- 
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FIG. 2: The dispersive behavior of the spin excitation in the superconducting state (8 — 0.125), in comparison with the 
spin- wave dispersion at half-filling (dashed curve). The peak positions of \" ln Q~ an d w-space are shown along different 
Q-directions: (a) along the diagonal direction, Q = (q,q); (b) along Q = (n,q). The solid bars mark the widths of the peaks 
in the momentum space (see text). 



Namely, there is no signature of x"(Q,c<j) at the ferromagnetic momentum Qo. 
At the AF momentum Qaf, one has 



C mm '(QAF) = ^X^ 1 ) 4 ^^VmatiVmvW^mvO')- 



(24) 



In the eigenequation (18), it can be easily shown that for any given state m there is a corresponding state m with the 
relation £ m = — £„ and uj ma {i) = (— l) 4 aj rf j Cr (i). Then (24) is reduced to 



Cw(QAF) = —^UJmai^lrAiWrn'A^rh'aii) = -fit 
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FIG. 3: Dynamic spin susceptibility x"(Qaf = (7r,7r),u>) in the superconducting phase with 8 — 0.125 (solid curve) 

the position of the resonancelike peak. The dotted curve is for the AF state at half filling. Inset: the evolution of the resonance 

peak at various dopings. 



and the dynamic spin susceptibility at Qaf can be simplified to 



x"(Qaf,u;) = ^ E (%) s § n H *(M - 2E ^- 



(25) 



The numerical result of x"(Qaf,^->) at 8 = 0.125 is shown in Fig. 3 by the solid curve. The dotted curve is 
calculated at half-filling, which diverges as 1/uj 2 at lo — > 0, in consistency with the spin-wave theory. Thus, in the 
superconducting phase, a resonancclikc peak appears at Qaf with a finite energy E g = 0.53J at 0.125 (twice bigger 
than that of E m shown in Fig.l). Note that higher energy (harmonic) peaks in x"(Qaf,w) are greatly reduced in 

strength in Fig. 3 due to the coherence factor -gj- in (25). So only the lowest peak at E g is clearly exhibited around 

Qaf- 

We further plot the resonancelike peak energy E g as a function of hole concentration in Fig. 4. At small doping, 
Eg is linearly proportional to 6, E g = 3.3(5 J, which is extrapolated to zero at half-filling, where the gapless spin wave 
is recovered. Note that in the present approach, the superconducting ground state is extrapolated to 5 = + . In a 
more careful study of the low-doping regime (beyond the mean-field approximation in the phase string model) has 
revealed that the AF state actually will survive up to a finite doping concentration, 6 < x c ~ 0.043 24 . In that case, 
one finds 24 that E g vanishes at 5 — x c following a square root behavir: E g ex y/S — x c as S — > x c , as shown by the 
dashed curve in Fig. 4. 

The momentum profile of the resonancelike peak at E g is shown in a three-dimensional plot in Fig. 5 at 5 = 0.125. 
It shows an intrinsic broadening of x"(Q,_E g ) in momentum around Qaf> which can be well fit by a Gaussian 
distribution function 



x"(Q,E g ) ocexp 



(Q 



Qaf) 2 



2a 2 



(26) 



The results for different hole concentrations are given in Fig. 6(a) along the diagonal momenta Q = (q,q). One can 
adjust a to make all data well collapse onto a single Gaussian function of (26) as shown in the inset of Fig. 6(a). 
The obtained broadening a turns out to be nicely scaled linearly withvo [see in Fig. 6(b)]. Similar plots can be done 
along different Q scans centered at Qaf and generally one has a oc vo in all directions. 




FIG. 4: The doping dependence of the resonancelike peak energy E g . The straight line illustrates the linear doping dependence 
at small 5. The dashed curve shows a E g oc y/S — x c behavior if the AF state survives at a finite doping x c , as shown in Ref. 24 . 
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FIG. 5: Momentum distribution of x"(Q,w) at u = E g (6 = 0.125). 



If we neglect the small anisotropy along different momentum directions centered at Qaf and perform a Fourier 
transformation to (26), we obtain the real-space correlation 



cr 2 R 2 R 2 

x"(R,w) oc exp( — ) = exp(--2-) 



(27) 
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FIG. 6: (a) Momentum distribution of x"(Q, E g ), scanned along the diagonal direction Q = (q, q) at various hole concentrations. 
The intensities are normalized at the maximums. The inset shows that the data in the main panel can be well fit into a Gaussian 
function exp(— (Q — Qaf) 2 /2ct 2 ), with a being scaled as linearly proportional tovJ, as shown in (b). 
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FIG. 7: The relation between the spin correlation length £/a and the hole concentration 8. The solid curve is —7=- The inset 
is the experimental results given in Ref. 9 . 
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with £ = *2.. Thus, the spin-spin correlation function decays exponentially with the distance in the superconducting 
phase. This is consistent with a spin gap E g opening up in the spin excitation spectrum. In Fig. 7, £ is well fit by 
the solid curve 



£ = «\/-7- 

7T0 



(28) 



In the inset, the experimental result obtained in LSCO 9 is presented for comparison. The general trend of £/a ex 1/\M 
in both the experiment and theory is quite telling. 

In Sec. Ill, we shall further discuss the momentum profile and longer spin correlation lengths at lower energies,, 
related to those seen in the LSCO compound 9 ' 10 , when the fluctuation effect is considered. 
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FIG. 8: Local susceptibility Xl( w ) i n the superconducting phase (solid curve) with S = 0.125, and at half filling (dashed curve) 
Inset: the evolution of the lowest peak (the resonance peak) at different dopings. 



4- Local susceptibility and spin sum rule 



The local spin susceptibility Xl( w ) 1S a ^ so an important quantity. It is defined by 



^>-/^«-«). 



(2*)' 

which describes the on-site spin-spin correlation. Based on (22), one obtains 



1 sgn(w)<J(|w| - E m - E m i), 



(29) 



where 



l*-mm' — ,,- / ^ yjJm<7\J j )\ \^m'a\ 



d\ 2 . 



(30) 



The numerical results of Xl( w ) at 5 = 0.125 and 5 = are presented in Fig. 8 by the solid and dashed curves, 
respectively. The low-energy parts in both cases are similar to those seen in \"(Q,uj) around Qaf (Fig. 3), as the 
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AF correlations are dominant at low energies. At high energies, more excitations which in momentum space disperse 
away from Qaf? as shown in Fig. 2, arc clearly present in x^(w). We see that the main band extends up to ~ 2.3J at 
half-filling, while is slightly reduced to around ~ 2.1J at 8 = 0.125. These upper-bound spin excitations are expected 
to be seen near the Brillouin zone boundary (see Fig. 2). 
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FIG. 9: The spin spectral weight vs doping. Open squares: the total weight; closed squares: the weight of the resonancelike 
peak around Qap; closed triangles: the weight of the AFLRO peak at half filling. The dashed line is the total weight from the 
exact sum rule, which is rescaled to coincide with the mean-field value at 8 = in order to compare the doping dependence. 



Although the intensity of each peak is physically not very meaningful, the weigh of the peak is. The reason is that 
there is a sum rule about the local dynamic spin susceptibility: 

Wtotai = J do, [l + n(w)]Xi(w) = ((S?f), (31) 

where the Bose distribution n(u) — 1/ (er u — l) . It means that the total weight of the spin susceptibility is related to 
an averaged spin number per site. At half filling, it is obviously that {(Sf) 2 ) is exactly 1/4. At finite doping, ((5f) 2 ) 
should be reduced (1 — 8) /A. 

In the bosonic RVB mean field state, the total weight can be calculated by 



W, 



total 



47V 4-^ 



1 + (blKr) 



(32) 



By using (b' ia bi a ) = (1 — 8)/2, we have W to ta] = |(1 — 8)(3 — 8). At half-filling, the total weight is | as compared 
to the exact result 1/4. The discrepancy is due to the relaxation of the no double occupancy to a global level in 
the Schwingcr-boson mean-field theory 20 . In Fig. 9, the doping dependence of Wtotai is shown with the exact result 
(dashed line) rescaled at 8 = 0. 
We also show the integrated weight of the resonancclike peak in Fig. 9 (solid curves with full squares), defined by 



Wpeak = 



peak 



du> [1 + n(u))j x'lM- 



(33) 



At 8 = 0.125, the weight of the peak is about 0.09, while the total weight is about 0.314, i.e., nearly 1/3 of the 
total weight is concentrated on the resonancclike peak. In Fig. 9, one can see that with the increase of doping 
concentration, W pca k actually gets slightly increased, whereas Wtotai is reduced. Namely, the resonancelike peak in 
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the superconducting phase will become even more prominent approaching the optimal doping from the underdoping. 
On the other hand, as the doping concentration is reduced to zero, IFpeak does not simply vanish. Instead, it 
approaches to a finite value which precisely coincides with the weight of the delta function at u> — and Q = Qaf 
in the dynamic spin structure function at half-filling, which represents the AFLRO. Earlier on, we have seen that at 
5^0 both Eg and the width a of the peak in momentum space go to zero. So the resonancelike peak continuously 
crosses over to the AFLRO at half-filling. 

III. SPIN DYNAMICS BEYOND MEAN-FIELD APPROXIMATION 

So far our discussions on spin dynamics have all been based on a generalized mean-field theory, characterized by 
the RVB order parameter A s [(12)]. Such a mean-field theory works quite well at half-filling over a wide range of 
temperature (~ J/ks) hi describing various ranges of spin-spin correlations. In particular, the nearest-neighbor (nn) 
spin correlation is directly related to A s by 



<Si-S i ) nn = --|A'|. (34) 



It thus provides an important justification for the doped case: Since spin-spin correlations, especially short-ranged 
ones, should not be 'washed out' immediately by the holes at small doping, the nn RVB pairing A s and thus the 
present mean-field state underpinned by the RVB order parameter is expected to persist over a finite range of doping, 
so long as the spin correlation length is no less than the nn distance (i.e., the lattice constant). In general, the effective 
Hamiltonian (10) is only valid within a low-doping regime of A s ^= which defines a pseudogap regime in the phase 
string model. Since a spin gap opens at finite doping in this regime, as shown in the last section, the amplitude 
fluctuation of the RVB parameter usually is not very important. 

Furthermore, we note that even within such a pseudogap phase characterized by a finite A s , the effective RVB 
description is not a usual mean-field theory beyond the half-filling. Generally speaking, the effective spinon Hamil- 
tonian (10) is a gauge model, in which the topological gauge field Afj describes w fluxoids bound to holes according 
to (11). Namely, this is not a spinon-only model and the hole-doping effect enters the Hamiltonian via A^, which 
represents the nontrivial frustration on the spin degrees of freedom from the motion of holes. In the previous section, 
the effect of A\, has been treated in a mean-field approximation. In the following, we shall discuss how to go beyond 
this mean-field level. 

A. Fluctuations induced by the density fluctuations of holes 

To examine the effect of fluctuations in A^ on spin dynamics beyond the mean-field approxiamtion, one has to first 
deal with the hole density fluctuations. In the phase string model, the hole degrees of freedom is also dependent 18 on 
the spin degrees of freedom. The nature of such mutually entangled charge and spin degrees of freedom is expected 
to make the theory quite nontrivial in a general case. 

In the superconducting phase, a uniform holon condensation 18 makes the topological gauge field A\, simplified as 
it may be treated as describing a uniform flux, namely, A\, ss A^ , with A^- determined by 

X;^- = *£*? = *£* (35) 

c lee lec 

for an arbitrary loop C according to (11). In the previous section, we have found that the spin dynamics in the 
superconducting phase is qualitatively modified by such A^ as compared to the AFLRO state at half-filling. 

However, the ideal Bose-Einstein condensation in treating A\, as A^- is only an approximate description of the 
holon condensation in the superconducting phase. In reality, one can expect all kinds of hole density fluctuations. 
The fluctuation of A\ ■, i.e., 6A^- = A\,— A^, will be tied to the density fluctation of the holes according to (11) as 
follows: 



E 



6AZ J =TrYM = irY(n?-6) (36) 



lec lec 



In the following, we shall examine the effect of 6A^- on the spin susceptibility previously obtained in the mean-field 
approximation. 
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Since holons are condensed in the superconduting phase, we may still reasonably neglect, to leading order of 
approximation, the dynamic fluctuations in SA^- and only focus on the static spatial fluctuations. As a Bose condensate 
is compressible, impurities and lattice distortions can all lead to some microscopic spatial inhomogeneity of the 
hole distribution, and below we introduce an approximate scheme to simulate 8A^- related to a microscopically 
inhomogeneous distribution of holes. 




FIG. 10: A way to introduce the gauge-field fluctuations related to the holon density. The open circle denotes a holon. The 
7r-fluxoid bound to the holon is smeared to the shadow area which is smaller than the whole lattice. 
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FIG. 11: x"(Qaf,^) with incorporating the fluctuations induced by the charge degrees of freedom. The doping is at 0.125, 
and the inset shows the local susceptibility Xi( w ) m the same situation. 



We first smear each 7r fluxoid bound to a hole within a finite size (Fig. 10 shows one configuration), representing 
some characteristic length scale of coherence for a bosonic holon, which should be still much larger than the average 
hole-hole distance to reflect the holon condensation. Then putting these smeared 7r fluxoids randomly on the lattice. If 
the smearing size of each 7r fluxoid is infinite, then the problem reduces back to the case of ideal Bose condensation with 



15 

5 A 1 ?, = 0. For finite sizes of fluxoids, there generally exist intrinsic fluctuations in the flux distribution of A\-, which 
we use to simulate the fluctuations related to the hole distribution. Since it is static, with each of such a configuration 
of non-uniform fluxes, we can follow the steps in last section to get a non-uniform mean-field solution and determine 
a dynamic spin susceptibility. The dynamic spin susceptibility at Qaf, averaged over the random configurations, is 
presented in Fig. 11, and the local susceptibility Xz,( w ) i s shown in the inset. The result is calculated in a 16 x 16 
lattice with each 7r flux being smeared within a 14 x 14 lattice size, with more than 10,000 configurations being 
averaged. 

For comparison, the mean- field results are plotted as dashed curves in Fig. 11. The main effect of such fluctuations 
in Afj is to cause the broadening of the resonancelike peak as well as high-energy peaks in energy space, although 
the peak positions, like E g , essentially do not change. Since in the mean-field case, the discrete levels are composed 
of degenerate Landau levels of £ m , a broadening due to lifting up the degeneracies by the fluctuations in A\- can 
be easily understood. So the above simple-minded approach to treat &A\, provides some valuable insight into the 
fluctuation issue in the framework of the bosonic RVB theory. A realistic treatment with a more accurate profile of 
the lineshape in the dynamic spin susceptibility is beyond the scope of this work. 

B. Incommensurability in momentum space 

The bosonic RVB mean-field state is based on the phase string formulation 19 of the t — J model, in which the 
short-distance singular part of the phase string effect introduced by the hopping of holes has been 'gauged away' such 
that the Hamiltonian in the new formalism is free of such singularities and thus becomes perturbatively treatable. 
But when one considers the physical quantities like the dynamic spin susceptibility, such singular effect should be 
still present in the correlation function and has to be incorporated carefully It has been shown previously 25 that 
the leading order contribution of such a singular effect to the dynamic spin susceptibility is simply represented by 
the incommensurate shifting of the momentum Q in C TO m'(Q) defined in (23) by 5Q X — ±2irg and SQ y — ±2irg 
(taking a = 1) with g ~ S. However, since the momentum width of the resonancelike peak in the mean-field is 
given by a — Vtt6, the incommensurability does not explicitly show up in the dynamic spin susceptibility 25 and the 
resonancelike peak still looks like one peak centered at Qaf, as illustrated in the top panel of Fig. 12. 

Now, due to the above-discussed fluctuational effect, the resonancelike peak is broadened with some of its weight 
shifting towards lower energies shown in Fig. 11. The corresponding width for these new low-lying modes in momentum 
space will be reduced too (i.e., the spin-spin correlation lengths are enhanced at energies lower than E g ) such that 
the incommensurability may become manifested in the dynamic spin susceptibility gradually with the decrease of the 
energy. Indeed, by using the same simulation used in Fig. 11, the incommensurate peaks do show up in the modified 
x"(Q>w) with incorporating the incommensurate shifting 25 , as u) is lowered below E g , which is illustrated in Fig. 12 
at S = 0.125. 

Therefore, the incommensurability in the dynamic spin susceptibility function is an intrinsic property of the phase 
string effect 25 . But its visibility crucially depends on spin fluctuations with longer correlation lengths at low energies. 
Such low-lying spin excitations, induced by the charge density fluctuations discussed above, are usually most prominent 
in the single-layer case, applicable to the LSCO compound. In contrast, the charge density fluctuations are expected 
to be weaker in the bilayer systems such as the YBCO compound, where the interlayer coupling will prefer the uniform 
distribution of the holons as to be discussed in the next section. 

IV. BOSONIC RVB DESCRIPTION WITH INTERLAYER COUPLING 

For a bilayer system, the t — J model can be generalized as 

^bilayer = _ f ^ ^.^ _ f± ^ C \ la Ci2a + H.C 
(ij)la i<? 

+J^S i; S j; + Jj_^Sa-S i2 (37) 

in which the additional subscript, I — 1,2, is the layer index. By introducing an additional bosonic RVB order 
parameter 



A s L =J2(b i i.b l2 -a), (38) 



and the Bogoliubov transformation 



biltr — ^^LJmka(i,l)( u mklmka — v mkl m k- a ) 



16 



(39) 



v i i A; 



Q=(7C+q. 7t) Q 




^ 




FIG. 12: The incommensurate structure is generally presented in x"(Q, w) due to the phase string effect 25 , but its visibility 
depends on the energy. The broad commensurate peak at u> = E g ~ 0.53 J at 5 = 0.125 is actually composed of four peaks 
(dashed curves), which becomes visibly 'split' in momentum space as u> is lowered below E g , when the fluctuational effect is 
included, where the individual peak width is reduced (or spin correlation length is enhanced) . 
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with k = ±, the mean- field spinon Hamiltonian can be diagonalized in the holon condensed phase in a procedure 
similar to Sec. II as given in Appendix A. We find 

Lo mk<J (i,l) = — = [sgn(£ m )fc] u ma (i), (40) 

(41) 

(42) 

V ^ V f^rak J 

and the spinon energy spectrum 




E mk = VA 2 - \Uk\ 2 , (43) 



in which 



Uk = sgn(U) [\U\ + -\^J ■ ( 44 ) 

In the above, £ m and to ma (i) are the solution of (18), as the counterparts of S, mk and u> mk(7 (i, /), respectively, in the 
single layer case. Finally, the self-consistent equations of the RVB order parameters and the Lagrangian multiplier A 
are given by 



|A S | 2 



VJZ^-R-T coth ^' (45) 



ANJ^ E mk 2 ' 



-(. k 



= ^ E g^ coth ^, (46) 

mfc 

2-5 = -A-^^coth^. (47) 

2N^-fE mk 2 v 7 

Fig. 13 shows the results obtained by the self-consistent equations as functions of doping concentration (solid lines) 
at J± = 0.1J, while the values of A s and A in the single layer case are plotted by the dashed lines for comparison, 
which only change slightly with the introduction of the interlayer coupling J± = 0.1J. 

In the inset of Fig. 13, the doping dependence of A^ is shown at various J±'s: J± = 0.10J, 0.11J, and 0.12J. 
We note that A^ is comparable with the interlayer pairing A s at half filling, e.g., 0.765 versus 1.157 even though 
J± = 0.1J is quite small. This may be attributed to the fact that the in-plane spin-spin correlation length £ is 
very large at half filling, which diverges at zero temperature. As the consequence, spin mismatches between the two 
layers will involve a large region determined by £, costing a big energy. This effectively enhances the interlayer AF 
correlations and thus the interlayer RVB pairing Aj_. 

Away from half filling, as shown in Sec. II, the in-plane spin-spin correlation length decreases monotonically with 
the hole concentration, which results in the reduction of the inflated interlayer AF correlations. Due to the competitive 
nature between A s and A^_ (one spin cannot be part of two RVB pairs at the same time), with the decrease of £, A^_ 
will diminish much faster than A s , as shown in the main panel of Fig. 13 as well as the inset for different Jj_'s. 

A. Spinon spectrum: Bonding and antibonding states 

According to (43) and (44), we find that with a finite A^_, the original spinon spectrum E m in the single-layer case 
is split into two branches, a bonding state E m+ and an antibonding state £" m -. 

The DOS of spinon spectrum at half-filling is shown in the inset of Fig. 14. At 6 = 0, £ m reduces to £ q = 

— JA s (cosq x a + cosq y a) . The ground state still has an AFLRO such that A = max(|£ q> fc|) = (2JA S -\ — =! ^- i )- The 
DOS of the bonding states is the same as the single-layer case at u) — > 0, while the antibonding states open a gap 
= min(_E m _) = 2^/JA s ^J±A s L ~ 0.60J at J± = 0.1J as shown by the dashed curve in the inset of Fig. 14. This gap 
is approximately the same as the gap in the dynamic spin susceptibility in the even channel (sec below) as observed 
by neutron scattering, which is about 70 mcV in magnitude 26 ' 27 . 
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FIG. 13: The doping dependence of Aj_, A s , A. The solid line is at J± = 0.1J, the dashed line is the results in the single layer 
case. The inset is the doping dependence of A s ± at different J±'s: solid line, J± = 0.10J; dashed line, J± — 0.11J, and dotted 
line, J± = 0.12J 

The main panel of Fig. 14 shows the DOS of spinon spectrum at 5 = 0.125 and J± — 0.1J in the superconducting 
phase. Compared to the spinon spectrum in the single-layer case with the Hofstadterlike structure illustrated by 
dotted lines in Fig. 14, there are bilayer splittings between the branches of the bonding (solid lines) and antibonding 
(dashed lines) states, given by E m - — E m+ . Because Jj_A^_ is much smaller than A, the splitting is most visible at 
the lowest energy level where £ m is the closest to A, as shown in the figure. In the following, we study how this bilayer 
splitting effect is manifested in the dynamic spin susceptibility. 



B. Dynamic spin susceptibility 



In the bilayer case, the imaginary part of the spin susceptibility \" [(Q, Q±) , w] depends not only the in-plane wave 
vector Q, but also the c-axis wave vector q±. It can be shown that 



x"KQ,q±)M =Xo(Q,w) sinV/2) + x' e '(Q,^)cos 2 (< Z± /2), 



(48) 



where x e * s ^ ne imaginary part of the spin susceptibility in the channels with odd and even symmetries, respectively, 
obtained from the retarded versions of the Matsubara Green's functions defined by 



Xo(i,r,r) = -(T 7 .(5f 1 (r)-5f 2 (T))(5| 1 (0)-^ 2 (0))) ) 
X e(i,j;r) = -(T r (5f 1 (T)+5f 2 (r))(^ 1 (0) + ^ 2 (0))). 



(49) 
(50) 



With the same procedure as in Sec. II C, we can obtain \o ano ^ Xe a ^ zero temperature straightforwardly as follows 



X "(Q,w) = ^ J2 C rnm>(Q)(l-Sgn(UU>)kk') 



>kk> 



Xe(Q: w ) = ^ X! C mm > (Q) (1 + sgn(£ m £ m ,)kk') 



32 



>kk' 



^mk^m'k' 

A <,mk^m'k' 

i^mkt^m'k' 



1 5(u> — E m k — E„ 



— 1 I 8(u> — E m k — E m 'k>). 



(51) 



(52) 
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FIG. 14: The DOS of spinons at 5 = 0.125 in the bilayer case. The solid line is for the bonding state, the dashed line for the 
antibonding state, and the dotted line denotes the single layer case for comparison. The inset shows the bilayer case at half 
filling: solid curve is for the odd channel and dashed curve is for the even channel. 



According to Sec. II C 3, C mm '(QAF) = jf^mm' such that 



SN 



ink 



mk 



x"(Qaf,w) 



-y 



167V 



^ ~r Qmkwm—k 
&'mk-&m—k 



i)6(w 



E m -k)- 



(53) 
(54) 



The above expressions clearly show that Xo(Qaf>w) is solely contributed by a pair of spinon excitations both from 
the bonding or antibonding states, while x"(Qaf,w) is contributed by a pair of spinon excitations, one from the 
bonding state and the other from the antibonding state. Compared to (25), one can see that Xo(Qaf) is very similar 
to x"(Qaf,u>) in the single-layer case. 

We present the numerical results at J± = 0.1 J and 5 = 0.125 in Fig. 15(a) and 6 = in the inset for comparison. 
The solid curve represents the odd mode while the dashed curve is for the even mode. From the main panel of 
Fig. 15(a), one sees that the single resonancelike peak in the single- layer case is replaced by a double-peak structure 
corresponding to the lowest bonding and antibonding states, respectively. In contrast, in the even channel, there is 
only one peak whose center is just in the middle of the double peaks in the odd channel. We also calculate the local 
spin susceptibility by integrating over the in-plane wave vector Q, which is given in Fig. 15(b). Two figures look 
quite similar. 

The doping dependences of the energies of these peaks are plotted in Fig. 16, where the filled squares mark the 
double peaks in the odd channel and the open squares describe the peak in the even channel. One finds that the 
doping dependences for the three peaks are very different. As S tends to zero, the lowest peak in the odd channel 
behaves like the resonancelike peak and reduces to the gapless spin wave mode at half filling, while the peak at a 
higher energy in the same channel move to high energy and reaches 1.190 J finally. 

At last, we consider the effect of the holon fluctuations in the bilayer case by using the same method introduced 
in Sec. Ill, and the results are plotted in Fig. 17, in which the solid curve is in the odd channel while the dashed 
curve is in the even channel. However, we point out that the interlayer coupling should be even more sensitive to the 
in-plane charge density fluctuations because the nonlocal phase factor involved in A^_ (see Eq.(A7) in Appendix A), 
which is not considered in Fig. 17. Generally speaking, the in-plane flux fluctuations due to the charge fluctuations 
will strongly frustrate the interlayer coupling. Thus, in the bilayer system, the former should be suppressed more, as 
compared to the single-layer case, by the interlayer coupling. More studies along this line will be conducted in future. 
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FIG. 15: (a) Dynamic spin susceptibility in the bilayer system at Qaf and S — 0.125. The solid line is for the odd channel and 
the dashed line is for the even channel. The inset is for the half filling case, (b) The local spin susceptibility. The solid curve 
is for the odd channel, and the dashed curve is for the even channel. 



V. CONCLUSIONS AND DISCUSSIONS 



In this paper, we have studied the spin dynamics in the superconducting state of a doped Mott insulator, which is 
described by the phase string model. In this model, the spin degrees of freedom are characterized by the bosonic RVB 
mean-field state, which can continuously evolve into the AFLRO state in the zero-doping limit, where the correct spin 
wave excitations are recovered. 

Our study has systematically shown how the low-lying spin wave excitations at half-filling are re-shaped into non- 
propagating modes in the superconducting phase by the motion of doped holes, via the phase string effect. We have 
found that the resonancclike peak near the AF wave vectors in the superconducting phase has its dominant spectral 
weight, at small doping, originated from that of the AFLRO at half-filling. That is, with the opening up a spin 
gap at finite doping, the low-lying spectral weight, including that of the condensed spinons, is pushed upward to a 
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FIG. 16: The energies of the peaks shown in Fig. 15(a) as functions of doping. Close squares: the double peaks in the odd 
channel; open square: the peak in the even channel. Inset: the difference between the peak in the even channel and the lower 
energy one in the odd channel. 



finite characteristic energy of the resonancelikc peak, which is linearly proportional to the doping concentration at 
small doping. We have analyzed the momentum broadening of the resonancelike peak, which decides a characteristic 
spin-spin correlation length, inversely proportional to the square root of doping concentration, or the average hole-hole 
distance. 

Our results have also clearly illustrated that the high-energy part of the dynamic spin susceptibility near the 
Brillouin zone boundary remains essentially the same at half-filling and at small doping, with the high-energy spin- 
wave signature still present in the superconducting phase. It reflects that fact that the local and high-energy AF 
correlations, within the length scale of the average hole-hole distance, have not been drastically changed by the 
motion of the holes. This is in sharp contrast to the prediction based on a Fermi-liquid- like theory, where the Fermi 
energy will serve as the natural high-energy cutoff in the spin susceptibility function. Our theory suggests that one 
must combine both the low-lying and high-energy spin excitations in order to correctly understand the nature of the 
spin dynamics in the high-T c cuprates. 

The fluctuational effects beyond the mean-field theory have also been examined. In the phase string model, the 
characteristic fluctuations will come from the density fluctuations of holons, which result in the local fluctuations of 
fluxes attached to holons but seen by spinons. The influence of such fluctuational effect on the spin degrees of freedom 
has been found to generally cause the broadening of the resonancelikc peak in energy space, making the emergence 
of some low-lying weight below the resonancelikc peak. This type of fluctuations is intrinsic and is believed to be 
important for the single-layer systems like the LSCO compound. In particular, we have found that the spin-spin 
correlation lengths of these low-lying modes are generally longer than the 'norm' one discussed in the mean-field 
theory. As a consequence, the incommensurability of the spin dynamics at low energies will show up, which is an 
intrinsic effect of the phase string model but is usually not visible when the width of each peak is too broad in 
momentum space, as in the 'norm' mode at the resonancelike peak. 

We have further investigated the interlayer effect on the spin dynamics by considering a bilaycr system. At half- 
filling, the spin excitation spectrum remains the same in the odd channel as the single-layer one at low energy, while 
a gap is opened up in the even channel, with the magnitude consistent with the experiment. Then we have shown a 
systematic evolution of the spin excitations, in both odd and even channels, with doping. In the superconducting phase, 
the effect of the interlayer coupling is most important for the low- lying resonancelike peak near Qaf- A prediction 
for the odd channel is that there will be a second peak with a smaller amplitude emerging at a higher energy, lying 
between the main resonancelikc peak in the odd channel and the peak in the even channel, near Qaf- However, both 
this second peak in the odd channel and the peak in the even channel will be sensitive to the fluctuations between the 
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FIG. 17: x",e(QAF,o>) with incorporating the in-plane holon fluctuations, simulated in the same way as in Fig. 11. The odd 
mode: solid curve; the even mode: dashed curve. 

two layers, which are not included in our mean-field treatment. 

Finally, we point out that in the present approach, our main efforts have been focused on the effective Hamiltonian 
H s , which describes the spinon degrees of freedom in the phase string model. The charge degrees of freedom are 
described by a holon Hamiltonian, Hh, in the phase string model, which is not considered explicitly as the holons are 
simply assumed to be Bose condensed in the superconducting phase, and two degrees of freedom are thus decoupled in 
this sense. But due to the mutual topological gauge fields in the phase string model, the condensed holons can feel an 
excitation from the spinon degrees of freedom and do response to it, as discussed in Ref. 28 . As a matter of fact, such a 
response will result in a loose confinement of spinons to allow only the S = integer types of spin excitations. We have 
considered the effect of Hh within the RPA and ladder-diagram approximations and found that the results presented 
in this work are not changed essentially, due to the fact that the interactions introduced by Hh are of logarithmic 
type and the spinons excitations are localized in space in the superconducting phase. Due to the length of the paper, 
we shall present these results in a separated paper. Lastly, we remark that the superconducting phase is not stable in 
the phase string model when the doping concentration is very low (< 0.04) where the spin ordered phase will persist, 
with the doped holes being localized, which have been discussed in Ref. 24 recently. 
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APPENDIX A: BOSONIC RVB MEAN-FIELD THEORY FOR THE BILAYER SYSTEM 

In the following, we generalize the bosonic RVB mean-field theory for the single-layer case 19 to a bilayer system as 
described by the generalized t — J model (37). 

We start with the phase string decomposition for the single-layer case 19 with explicitly introducing a layer index I 
for each layer (I = 1,2): 



j. -.--.string 



(Al) 



2: J , 



where e <■'" tracks the in-plane phase string effect, defined by 



&2r* = l{* b u-<r*b, (A2) 



with 









a»Sia - 1 . ( A3 ) 



and 



The exchange term in the phase string formulation reads 

F bi. ayer = _ J ^ (A^) f A^ - ^ ]T (A| i;J .) f A& fX . (A5) 

where the in-plane RVB pair order parameter 

A ?i,i=E e "" 4 ' ,!) ^ 6 i'- ( A6 ) 

and the intcrlaycr RVB pair order parameter 

^,±=E e " <f( *""* f!,)6 ^^-«r. (A7) 

In the single-layer case, the hopping term contributes to an additional feedback effect 18 on the spin degrees of freedom, 
besides the phase string effect. But it does not qualitatively and quantitatively change the main results of the spin 
dynamics. Similarly, the interlayer hopping term is not considered here due to the same reason. 

In the superconducting state, due to the holon condensation 18 , the in- plane gauge field A\- 1 can be treated as 
describing a uniform flux [cf. (14)]. On the other hand, the phase difference between 3?^ and $^ 2 for two layers may 
be considered as a constant in the holon condensation case, i.e. $^i — $^2 = 4>i so that it can be gauged away. Then 
it is natural to introduce the following RVB order parameters 

A s = (^e-^^'biurba-a), (A8) 

Ai = (^bnvba-a). (A9) 

The superexchange term including the interlayer coupling is thus reduced to 

(ij)cjl ia 

+A(E 6 L^-2(l-*)iVj. (A10) 

\ila ) 

To diagonalizc this Hamiltonian, we introduce the generalized Bogoliubov transformation 



bilcr = ^ j U m k<j{i,l){u m klmk<7 ~ VmklJnk-a), (AH) 

where k = ±. By requiring 



i k 



t 1 _ _ TP . J 



[H s ,Jmka] = E mklmka , and [ff s ,7„J = -£ ra fc7™to> ( A12 ) 
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wc find 



(\-E mk )u mkl7 (i,l) = -^ E ^-^O^-^^^-^O (A13) 



j=NN(i) 

s 

2 " 



(A + £U)tw(M) - ~ E ^fc-^O'.Oc^^-^y^^fc-^i.n. (A14) 

j=NN(i) 

where I' denotes the layer different from I. We obtain the solution 

u mktT (i, I) = U mk LO mka {i, I), v m ka{i, I) = v mk aj mklJ (i, I), (A15) 

with 

«L - "L = !. ( A16 ) 

and u mka {i,l) satisfies 

imk^mU%l) = -^ E ^mk-AJJ)e i(TA ^-^Y L ^ mk -Ai,n- (A17) 



2 

j=ATJV(i) 



The spinon spectrum is given by 



and 



£„,» = , a 2 -e** (Ais) 



"- = vK^i +1 i 



According to (18), one has 



v mk = sgn(£ mk )J- I— f). (A20) 



Vmka(i,l) = —7= (Sgn(£ m )fc) U>mer(i), ( A2 1) 

g mfc =sgn(^ m ) |g m | + x , (A22) 



where £ m and o; mcr (i) are the solutions of (18). 

Finally, the self-consistent equations of the RVB order parameters and the Lagrangian multiplier A can be obtained 
as follows 



|A'T = — y^^coth^^, (A23) 



2 



4AJ^f E mk 2 



Sgn(£ m )fc£ m fc ,, PE, 

cotli 

2-5 = ^y^coth^^. (A25) 



AI = _ 1 ygg%lM^ coth ^^, (A24) 

± 2N ^ E mk 2 V ; 



27V ^— ' i£ m /g 
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